import numpy as np
import matplotlib.pyplot as plt
import math
import pandas as pd
import statsmodels.api as sm
from scipy.stats import ortho_group
import random
np.random.seed(7)
def f(x1):
    return sum([-label[j]*np.dot(A[j],x1)+math.log(1+math.exp(np.dot(A[j],x1))) for j in range(la)])
def bijing(x1,delta,n):
    v = ortho_group.rvs(dim=ss)
    bb= np.zeros([ss,1])
    for j in range(ss):
      vb = v[:, j]
      vb = vb.reshape((ss, 1))
      cc = (1/(2*delta))*(f(x1+np.dot(delta,vb))-f(x1-np.dot(delta,vb)))
      bb = bb + np.dot(cc[0], vb)
    return bb
def grad(x1):
    D = [sum([-label[j]*A[j,i]+A[j,i]/(1+math.exp(-np.dot(A[j],x1))) for j in range(la)]) for i in range(ss)]
    D = np.array(D)
    return np.reshape(D,[ss,1])
def NAG(origin, rate, p, epochs):
    x = np.array(origin)
    y = np.array(origin)
    z = np.array(origin)
    dot_set = [x]
    for i in range(epochs - 1):
        z = y
        y = x - np.dot(rate, grad(x))
        x = y + (1 - co(p, i)) * (y - z)
        dot_set.append(x)
    return dot_set
def co(p,t):
    return p/(t+p)
def compare(x1,x2):
    return x1 if f(x1)<=f(x2) else x2
#比较函数，输入两个n维list,输出一个n维list
def GNSA(origin,epochs,rate,p):
    x = np.array(origin)
    y = np.array(origin)
    z = np.array(origin)
    n = len(x)
    dot_set = [x]
    delta = 1
    for i in range(epochs - 1):
        if i <50:
           delta = delta*(0.5)
        y = (1 - co(p, i)) * x + co(p, i) * z
        x = compare(y - 2*rate * bijing(y,delta,n),x-2*rate*bijing(x,delta,n))
        z = z - (rate / (co(p, i))) * bijing(y,delta,n)
        dot_set.append(x)
    return dot_set
def pr_fun_value(y,min):
    plt.semilogy([i for i in range(len(y))],[f(y[i])-min for i in range(len(y))])
    return
def GD(origin,epochs,rate):
    x = np.array(origin)
    dot_set = [x]
    for i in range(epochs - 1):
        x = x - rate * grad(x)
        dot_set.append(x)
    return dot_set
def QGNSA(origin,epochs,rate,p):
    x = np.array(origin)
    y = np.array(origin)
    z = np.array(origin)
    n = len(x)
    dot_set = [x]
    delta = 1
    for i in range(epochs - 1):
        delta = delta *((0.5)**(0.1))
        y = (1 - co(p, i)) * x + co(p, i) * z
        x = compare(y - 2 * rate * bijing(y, delta, n), x - 2 * rate * bijing(x, delta, n))
        z = z - (rate / (co(p, i))) * bijing(y, delta, n)
        dot_set.append(x)
    return dot_set
mm = 300
cend = 600
oo = 0.0002
pp = 3
la= 200
ss= 5
a = np.ones((ss,1))
a = np.ones(ss)
a = a.reshape(ss, 1)
A = np.random.normal(size=(la, ss))
label = np.reshape(np.array([random.randint(0,1) for t in range(la)]),[la,1])
cp = NAG(origin=a,epochs=cend,rate=oo,p=4)
min1 = np.min([f(cp[i]) for i in range(cend) ])
pr_fun_value(GD(origin=a,epochs=mm,rate=oo),min=min1)
pr_fun_value(NAG(origin=a,epochs=mm,rate=oo,p=3),min=min1)
pr_fun_value(NAG(origin=a,epochs=mm,rate=oo,p=4),min=min1)
pr_fun_value(QGNSA(origin=a,epochs=mm,rate=oo,p=3),min=min1)
pr_fun_value(QGNSA(origin=a,epochs=mm,rate=oo,p=4),min=min1)
#pr_fun_value(GNSA(origin=a,epochs=mm,rate=oo,p=3),min=min1)
#pr_fun_value(GNSA(origin=a,epochs=mm,rate=oo,p=4),min=min1)
plt.legend(loc='best',frameon=True, labels=['GD','NAG(p=3)','NAG(p=4)','NSA-zero(p=3)','NSA-zero(p=4)'],fontsize='large')
plt.xlabel('iteration')
plt.ylabel('f - f*')
plt.show()





























